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In the density functional (DF) theory of Kohn and Sham, the kinetic energy of the ground state 
of a system of noninteracting electrons in a general external field is calculated using a set of orbitals. 
Orbital free methods attempt to calculate this directly from the electron density by approximating 
the universal but unknown kinetic energy density functional. However simple local approximations 
are inaccurate and it has proved very difficult to devise generally accurate nonlocal approximations. 
We focus instead on the kinetic potential, the functional derivative of the kinetic energy DF, which 
appears in the Euler equation for the electron density. We argue that the kinetic potential is more 
local and more amenable to simple physically motivated approximations in many relevant cases, and 
describe two pathways by which the value of the kinetic energy can be efficiently calculated. We 
propose two nonlocal orbital free kinetic potentials that reduce to known exact forms for both slowly 
varying and rapidly varying perturbations and also reproduce exact results for the linear response of 
the density of the homogeneous system to small perturbations. A simple and systematic approach 
for generating accurate and weak ab initio local pseudopotentials which produce a smooth slowly 
varying valence component of the electron density is proposed for use in orbital free DF calculations 
of molecules and solids. The use of these local pseudopotentials further minimizes the possible errors 
from the kinetic potentials. Our theory yields results for the total energies and ionization energies 
of atoms, and for the shell structure in the atomic radial density profiles that are in very good 
agreement with calculations using the full Kohn-Sham theory. 
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I. INTRODUCTION 

Density-functional theory (DFT) has become one of 
the most powerful tools for investigating the electronic 
structure of large complex systems. In principle, as 
shown by Hohenberg and Kohn [ij], the exact ground 
state energy of a system of N electrons can be formally 
written as a functional E[p] of only the electron density 
p(r), a function of three variables, and the external field 
V ext {r). Determining the energy and other ground state 
properties from such an approach could dramatically 
reduce the computational cost for large systems when 
compared with traditional quantum chemistry methods, 
which deal with wavefunctions involving coordinates of 
all N electrons HQ. 

Kohn and Sham (KS) |, Q showed that E[p] can be 
usefully partitioned into the following set of terms: 

E[p] = T s [p] + E H [p] + E xc [p] + [ p(r)V ext (r)dr. (1) 



equals p(r) 



Here T s [p] is the noninteracting kinetic energy density 
functional (KEDF), which gives the kinetic energy of 
a model system of N noninteracting electrons in a self- 
consistent field chosen so that the ground state density 
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Eh[ P ] = - 



p(r)p(r') 



drdr' 



(2) 



is the classical electron- electron potential energy (Hartrce 
energy) and E xc [p] is the exchange- correlation energy (in- 
cluding the difference between the interacting and non- 
interacting kinetic energy and the difference between 
the quantum and classical electron-electron potential en- 
ergy) . The last term on the right of Eq. (fTJ) is the only 
term that depends explicitly on the external potential 
V ex t(r). Atomic units are used throughout the paper. 

If all these functionals were known, then the density 
p(r) could be obtained from the variational principle (Eu- 
ler equation) associated with minimizing Eq. (fTJ) : 



p = V Ts (r; [p]) + V eff (v; \p]), 



(3) 



and the total energy of the inhomogeneous system could 
then be determined from the energy functional E[p\. All 
other physical quantities related to the ground-state den- 
sity could also be computed. Here p is the chemical po- 
tential (the Lagrange multiplier associated with the nor- 
malization condition J p(r)dr — N), and V e ff(r; [p]) is 
an effective one-body potential defined by 



Veff(r; [p]) 



Sp(r) 
V H (r; [p]) 



E H [p]+E xc [p]+ J p(r)V ext (r)dr 
V xc (r- [p])+V ext (r), 



(4) 



where 



V H (r;\ 



6E H [p]/Sp(r) 



P(r') 



:dr' 



(5) 
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is the Hartree potential, and V xc (r; [p]) = SE xc [p]/Sp(r) is 
the exchange-correlation potential. Similarly we interpret 

V Ts (r;[p}) = ST s [p]/5p(r) (6) 

as the kinetic potential (KP) arising from the KEDF 6]. 

Further progress requires an accurate determination 
of the noninteracting kinetic energy, whose magnitude is 
much larger than the exchange-correlation energy. The 
initial development of DFT as a practical computational 
method was made possible by KS's realization that the 
numerical value of the noninteracting kinetic energy can 
be exactly calculated, not directly from the density itself 
using T^p], but by introducing a set of N one-electron 
wave functions (orbitals) satisfying the N coupled KS 
equations that describe the model system 0, || . 

Research could then focus on determining the remain- 
ing small term E xc [p]. Here even local density approx- 
imations have often proved useful. Through the efforts 
of many workers we now have generally accurate expres- 
sions for E xc [p\. Using these along with the KS orbitals 
to calculate the kinetic energy, one can accurately cal- 
culate both the total energy E[p] and the ground-state 
density p(r) for a wide variety of systems. 

However, the use of the KS orbitals usually generates 
a relatively expensive 0(N 3 ) scaling of computational 
cost with the number of electrons. While this scaling is 
much better than that of most standard methods that 
include correlation energy calculations for large systems 
remain problematic. This remaining bottleneck could be 
removed if there were an accurate treatment of the kinetic 
energy in terms of the electron density only 0, 0, 0] . 

To that end there has been considerable effort in- 
vested in developing "orbital-free" density functional the- 
ory (OF-DFTVby making direct approximations for T s [p] 

B 0, M HO, lLI El M M, M \m M, M Hi while 

earlier simple local density approximations for T s [p] like 
those used in the Thomas-Fermi (TF) model [2l[ are very 
inaccurate, there have been two main advances in recent 
work that offer prospects for significant improvements. 

The first is the introduction of nonlocal KEDFs that 
reproduce known exact results for very slowly varying or 
very rapidly varying fields and give the exact linear re- 
sponse (LR) of the density of the uniform model system 
to small perturbations. Similar ideas have been success- 
fully applied to classical nonuniform fluids [12] • The sec- 
ond advance is to focus not on the total density but on 
the smaller and more slowly varying valence electron den- 
sity as described by a weak pseudopotential acting only 
on the valence electrons. While conventional pseudopo- 
tential methods use orbitals, recently developed ab ini- 
tio local pseudopotential (AILPS) methods determine the 
unique local one-body potential producing a given target 
valence density by solving the KS equations inversely, us- 
ing the one-to-one mapping between density and poten- 
tials in DFT H|. For OF-DFT with LR-based KEDF's, 
the use of pseudopotentials not only can reduce the com- 
putational cost, but also can improve its accuracy, since 



the system will be closer to the LR regime where T s [p] is 
designed to be accurate. 0, . Indeed, very promising 
results using such OF-DFT methods have been obtained 
for a variety of nearly free-electron-like metals. 

However, existing KEDF's have not yet achieved chem- 
ical accuracy for systems with localized and more rapidly 
varying electron densities like molecules or for covalent 
or ionic solids. The main problem is that the exact T s [p] 
is highly nonlocal, and we have little idea of the func- 
tional form of the nonlocality for densities far from the 
LR regime. It has proved very difficult to understand 
what errors an approximate nonlocal T s [p] will produce 
in the density as determined by the Euler equation with 
a general V ext (r). 

We explore here a different way to attack this basic 
problem. The exact T s [p] can be formally obtained from 
Vt b (r; [p]) by functional integration over density changes 
in all regions of space [f| [H| . Because of this integration 
T s [p] is a more nonlocal functional of the density than 
is Vt, (r; [p]). More detailed arguments arriving at this 
same conclusion have been recently presented [26[ . Since 
most problems in devising accurate approximations for 
T s [p) have arisen from the nonlocality, this suggests it 
could be worthwhile to try to develop approximations 
for the KP V Ta (r; [p]) itself. 

To illustrate this point, Chai and Weeks Q added a 
simple gradient correction to the original local TF KP 
for atoms [2lJ, with a coefficient chosen to reproduce 
the exact boundary condition of exponential decay of the 
electron density far from the nucleus. Though quanti- 
tative results were not obtained, the resulting modified 
Thomas- Fermi (MTF) model gave energies for atoms and 
for closed-shell diatomic molecules that showed notable 
improvements when compared to the original TF and re- 
lated gradient corrected KEDF models. However, the 
local gradient correction used in the MTF KP cannot re- 
produce the oscillatory atomic shell structure and it does 
not satisfy the exact LR behavior in the homogeneous 
limit. It is clear that nonlocality even in the KP must be 
taken into account to achieve more accurate results. 

We propose here new nonlocal approximations for the 
KP using ideas similar to those employed for the nonlocal 
KEDFs. These new KPs satisfy the exact LR condition 
in the uniform limit, and reproduce known exact limiting 
forms of Vt s (r; [p]) both for very slowly varying and very 
rapidly varying perturbations. As will become clear, the 
nonlocality in our KP is determined by the requirement 
that LR is exactly satisfied, and it is much easier to en- 
sure that LR holds for the KP than it is for analogous 
KEDF models. We believe this level of nonlocality in the 
KP may suffice in many cases when used in conjunction 
with AILPS methods to describe slowly varying valence 
density components closer to the LR regime. 

The remainder of this paper is organized as follow. 
Section II will discuss some general pathways connect- 
ing T s [p] and Vr s (r; [p]). Section III will describe limit- 
ing forms of the KEDF and KP for slowly varying and 
rapidly varying perturbations, and discuss LR theory, an 
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exact theory for the response of density of the uniform 
electron gas to small perturbations. Section IV will de- 
velop two nonlocal KPs incorporating both the correct 
limiting forms of the exact KP and the exact LR of the 
free-electron gas. Section V will compare the numerical 
results of the present method for atoms with the KS-DFT 
and other KEDFs, both for all-electron calculations and 
for valence electrons using the AILPS. We find that the 
use of AILPS indeed reduces errors arising from nonlocal- 
ity in these approximate KPs or KEDFs, which give very 
accurate results for the relatively slowly varying valence 
densities. Our conclusions are given in Section VI. 

II. PATHWAYS FROM V Ts (r; [p]) TO T s [p] 

If T s [p] is known, Vr s (r; [p]) can be simply computed 
by functional differentiation. However, there is no unique 
way of determining T s [p] from a given Vr s (r; [/>])■ Many 
possible pathways can be used to construct T s [p] by func- 
tional integration of Vr s (r;[p]) [2?], [H)]. If the exact 
Vr s (r; [p]) is used and the integration is carried out ex- 
actly, then all pathways would give the same exact result 
for T s [p\. However, when an approximate Vr s (r;[p]) is 
used, different pathways will give different results for the 
kinetic energy. But this "thermodynamic inconsistency" 
is small if reasonably good approximations are used, since 
the integration tends to smooth out local errors that may 
exist in the density 0, [H[ • 

More problematic is the fact that most pathways re- 
quire additional results for partially coupled systems as 
the external field or density perturbation is gradually 
turned on, which adds to the computational burden. In 
particular, most earlier work has used a "potential en- 
ergy pathway", where the external potential is scaled by 
a coupling parameter The kinetic energy can then 

be found by subtracting the potential energy (calculated 
from the potential energy density functionals) from the 
total energy. However, this pathway is expensive, since 
one has to solve the Euler equation ^ for each partially 
coupled V* xt (r) (with the same p), to determine the cor- 
responding p\. 

A. Herring's Pathway 

However, Herring showed there is a particular path- 
way arising from exact scaling relations between the non- 
interacting kinetic energy T s [p] with respect to the coor- 
dinate r in p(r) where very simple results involving only 
the final density can sometimes be found [f| [l8| . If the 
coordinate r is scaled to ar, the normalized scaled den- 
sity is /9 Q (r) = a 3 p(ar). It is easy to show that the exact 
T s [p] then obeys 

a 2 T s [p] =T s {p a }. (7) 

For isolated systems, such as atoms and molecules, the 
density and its derivatives to all order vanish far from the 



nuclei. For such systems, when Eq. ([7]) is differentiated 
with respect to a, and the partial derivative is evaluated 
at a = 1, we find the formally exact result 

T M = \ J V Ts (r;[p])V.(rp(r))dr. (8) 

Therefore, once the kinetic potential Vr s (r;[p]) is 
known for some given p(r), the numerical value of T s [p] 
can then be immediately determined from Eq. ([5]) . Since 
there is no need to perform a coupling parameter integra- 
tion over the change of density or potential, this scheme 
is not only fast, but also numerically reliable. The final 
form of Eq. ([5]) is essentially the viral theorem, and is 
directly related to the force on molecules [29L l30l. l3ll. |32| . 

Note that this simple and exact pathway holds only 
for the noninteracting T s [p] [3~H |. which again shows the 
virtues of the KS partitioning of the total energy. We will 
use Eq. ((5|) as the basic pathway to determine the numer- 
ical value of T s [p] from a given approximate Vr s (r; [p]) 
for most calculations in this paper. However, Eq. ([5]) 
does not hold for extended solid-state systems because of 
nonvanishing boundary terms, and thus far we have not 
found an exact and simple way of including them. 

Fortunately, there is another class of computationally 
efficient "density pathways" that can be used for extended 
systems, as we now show. Density pathways can also 
be used for atomic and molecular systems to check the 
accuracy of the Vr s used, since results using the exact 
Vt s would be independent of path [2?], HH . 

B. Density Pathways 

The change in the kinetic energy can be formally re- 
lated to a coupling parameter integration, where the den- 
sity changes from some known value at A = to the final 
density at A = f : 

T s [p] = T X=0 +f Q d\j drV Ts (r; [px])^^- (9) 

In most cases a simple linear density pathway will suffice. 
Here the density p(r) is linearly scaled by a coupling pa- 
rameter A from some uniform reference density po natu- 
rally chosen to be the uniform electron density N/V in 
extended systems: 

px(r) = p + X[p(r)-p ]. (10) 

Then Eq. ((9]) becomes 

T s [p] = T x=0 + J dxj dvV Ts (r; [px])W) - Po] (11) 

Here T\=o is the kinetic energy of the uniform system, 
i.e., the Thomas-Fermi kinetic energy Ttf[po\- F° r ex ~ 
tended systems, where the Herring's pathway cannot be 
used, this density pathway appears to be a good way to 
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compute T. Other density pathways, like the square-root 
pathway introduced by Chen and Weeks [28[ to describe 
nonuniform hard sphere fluids, can be defined, and have 
proved useful in certain applications, but we do not con- 
sider them here. 

Note from Eq. (fT0)l that p\(r) depends only on the fi- 
nal density, so evaluation of Eq. (jlip is straightforward 
and this pathway is computationally efficient. Unlike the 
potential energy pathway where the external potential is 
scaled, there is no need to solve the Euler equation © 
for its corresponding external potential V* xt (r) at each A. 
However, for isolated systems, where po — 0, this path- 
way is likely to be less accurate than Herring's pathway, 
since it does not automatically satisfy the virial theorem. 



III. EXACT LIMITS AND LINEAR RESPONSE 
THEORY 

Although the exact T s [p] is still unknown, several limit- 
ing forms have been discovered for particular density dis- 
tributions. These provide important cornerstones that 
can be used to construct accurate KEDFs and KPs in 
many cases, as will be seen below. 

In particular, the Thomas-Fermi (TF) KEDF HJ is 
known to be exact for a uniform system: 



T TF [p] = C F I P b ' z {v)dv 



where C F = ^(3tt 2 ) 2 / 3 . The TF KEDF T TF [p] is de- 
rived by local use of uniform free-electron gas model, and 
is exact for a system with an infinite number of electrons. 
The corresponding expression for the TF KP is 

V TF {v- [p]) ee 5T TF [p]/5p{v) = \c F p 2 '\v). (13) 

This depends only on the local value of /3 2//3 (r) and thus 
formally is more local than the TF KEDF, whose func- 
tional dependence on p involves the density at all r. Of 
course in this simple case the functional integration of 
Eq. (fl3]) can be carried out exactly to yield Eq. (p~2|) . 
but this cannot be done in general and the nonlocality of 
T s [p] has proved problematic. 

Results for nonuniform systems are best described in 
Fourier space. For a very slowly varying perturbation of 
the density, the second-order gradient expansion is exact 
[33l ] . It is easy to see that results correct to second order 
at small wavevectors are given by 



where 



t tf±w[p\ = Ttf[p\ + gT w [p], 



p(r) 



is the von Weizs&cker (W) KEDF [34 



(14) 



(15) 



7V[p] is exact for a system with one or two electrons, 
or where the density can be accurately described by a sin- 
gle orbital. Moreover it has been argued I,©, El El that 
^V[p] gives the correct leading order term for a rapidly 
varying perturbation with only high wavevector compo- 
nents and that the next order correction is reproduced 
by 



T w _3 TF [p\ = T w [p] - -T TF [p\. 



(16) 



The W KP is 



Mr;W ) sJ r„ H / Wr )4(J™W!!- 2 ^M 

(17) 

If we represent the full density by an effective single or- 
bital function ip(r), 



p(r) = |V(r)| 2 



(18) 



then the W KP can be written in a compact form that 
will later prove useful: 



Vw(r; [p]) 



vV(r) 
2^(r) • 



(19) 



Finally, the linear response of the density of a uniform 
non-interacting electron gas with density po to a small 



(12) perturbation SV(k) = eke is exactly known [351 ] . 



Here 



Sp(h)= XL (q)SV(k). 



q ee k/2k F 



(20) 



(21) 



is a dimensionless wavevector, where 

k F EE (3^ 2 P0 ) 1/3 

is the Fermi wavevector (FWV) and k = |k| 
function xl(q) has the form 

XL (q) = -%F^{q) 



(22) 
The LR 



where 



FUq) 



7T 

hp 

'7? 



1 - q 2 



In 



1-3 



■In 



1-9 



-1 -1 



(23) 



(24) 



has been called the Lindhard function [10 ]. 

It is known that the weak logarithmic singularity at 
q = 1 in F^ L 1 {q) is responsible for Friedel oscillations, and 
may also be important for the appearance of atomic shell 
structure. This singularity further divides the Lindhard 
function into two branches in Fourier space: the low- 
momentum (q < 1 ) or the low-g (LQ) branch, and the 
high-momentum (q > 1) or the high-g (HQ) branch [lol ]. 
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Figure 1: Linear response functions of a uniform system of 
noninteracting Fermions as given by the TF, W, and MTF 
(see Ref. 0) models. 




The dimensionless response function arising from the 
TF KEDF is PW(g) = 1, and that from the W KEDF is 
F\V (?) = 3<7 2 [36[ . Clearly, no linear combination of the 
TF and the W KEDFs can reproduce the exact Lindhard 
function in Eq. ([24]) . This has the following two limits 



10], 



Fl W-\ 3g 2 -3/5 



0(q 4 ) 
^0(q- 2 ) 



q > 1 



(25) 



It should be noted that the expansions for both the 
low-q and high-q limits are correct to all orders in per- 
turbation theory, but valid only in the appropriate limits 
in Fourier space. On the other hand, the LR theory is 
valid for all wavevectors, but is only accurate for small 
perturbations. Therefore, the regime where the response 
functions of the two limiting KEDFs deviate from the 
exact LR function gives an indication of the range of 
wavevectors where the two limiting forms are inaccurate. 

As shown in Fig. [TJ the response function F^p{q) has 
no momentum dependence and is only exact at q = 0. 
The response function F^r (q) is exact asymptotically at 
high q, and remains fairly accurate for q > 2, but is di- 
vergent in the low-q branch, and fails completely for the 
nearly uniform electron gas. In contrast, the MTF model 
Q gives a reasonably accurate average description of the 
exact response function, especially in the important re- 
gion near the singularity at q = 1. 



IV. CONSTRUCTION OF NONLOCAL 
KINETIC POTENTIALS 

A. Kinetic Energy Density Functionals Ttf\w[p] 
and T w> ,tfIp] 

Simple linear combinations of the two limiting KEDF's 
in Eqs. (JUD and (TH), such as the TFAW KEDF 33, 37, 



Ttf\w[p] = T TF [p] + \T w [p] 
and the WATF KEDF gO, EH, E3, [H 

Tw\tf[p] = T w [p] + XT TF [p] 



(26) 



(27) 



have been widely studied for several decades. The value 
of the parameter A was either determined empirically for 
getting good atomic energy or obtained by some semi- 
classical arguments. 

The advantage of these approaches is the ability to gen- 
erate a family of simple KEDF's easily. It has been shown 
empirically that the TF1/5W model can give good val- 
ues for atomic energies, but the predicted density profiles 
are generally not very accurate, both near and far away 
from the nucleus. The Ttf\w[p] and Tw\tf[p] function- 
als give the correct leading term in the density response 
to a slowly-varying perturbation and a rapidly-varying 
perturbation respectively, and with particular choices of 
A as in Eqs. ([14]) and ([TO]) they can reproduce the next 
order term. Unfortunately, they then will have an in- 
correct leading term in the opposite limit, unless A = 1. 
However it has been shown that Ttfw[p] with A = 1 
always overestimates the exact T s [p] for various systems 
[nj. Finally, none of these functionals can reproduce the 
exact response function Fl(c[) in the homogeneous limit. 
Since these models fail to satisfy all the known limits, 
and nonlocality in T s [p] is not correctly described, it is 
also not surprising that atomic shell structure is missing 
in these approaches. 



B. Combining TF and W Kinetic Potentials 

We argue that it may be more profitable to take advan- 
tage of known limiting forms of the KP, rather than the 
KEDF, and develop approximations for the more local 
Vr s (r; [p]) directly. Again we can rely on known results 
in the linear response regime when the density variations 
are not too large. 

From Eqs. ([14]) and p6[) . the following linear combina- 
tions of the TF KP and the W KP in Eqs. ([13]) and (17]) 
can reproduce exact results to second order for very small 
and very large wavevector perturbations respectively: 



v T .(r;\p]) 



V TF (r; [p]) 
Vw(r; [p]) - 



\V w {v- [p]) 
\Vtf{t;\p]) 



q < 1 
q > 1 



(28) 



Since Vtf{y] [p]) and VV(r; [p]) are the only components 
up to second order of the two exact limiting forms of the 
KP, we can combine them in analogy to the TFXW and 
WXTF models and arrive at generalized KPs. 

However, instead of combining them using a fixed pa- 
rameter A, it seems natural to represent them in Fourier 
space and allow a wavevector dependence in A — X(q) to 
connect the limiting forms. The X(q) can then be chosen 
in a very simple way so that the exact LR function is 



6 



Figure 2: Weight function f(q) for the HQ and LQ KPs. 



reproduced for a uniform system with density po . In this 
way the LR function bridges the exact limits at large and 
small wavevectors, and if the theory is applied to weak 
perturbations in the linear response regime for interme- 
diate wavevectors we can expect very accurate results. 
Here, we derive such generalized KPs based on the KP 
for the WXTF model. 



C. HQ Kinetic Potential 

In analogy to the WATF model in Eq. ([2"T]) we look for 
a kinetic potential of the form 



where 



W(k) + Xhq(q) Vtf (k) 
VtfQ*) + VwQl) + f(q) Vtf (k) 

q = k/2k F 



(29) 



(30) 

is a dimensionless wavevector normalized by the FWV k F 
in Eq. (|22[) of a uniform reference system with density 
Pq and f(q) = \HQ(q) — 1- The superscript in V^q 
indicates use of a uniform reference system. For a small 
perturbation, we can linearize the V^g(k) in Eq. ([2"9"]) . 
Requiring that it satisfy LR exactly then determines the 
weight function f(q) as 



f(q) =F L (q)-3q 2 



1. 



(31) 



See Fig. H 

We refer to Eq. ([29]) with Eq. ([31]) as the HQ KP model. 
It reproduces the correct high-g limit in Eq. ([2"5]) up to the 
second order. However, unlike Eq. (fTB)) . it also satisfies 
the correct low-g limit to leading order and gives exact 
results for all q in the linear response regime. Inverse 
Fourier transform of Eq. ([2"9"]) then gives 



V° Q (r; [p], k F ) = V TF (r; [p]) + V w (r; [p]) 
+J /(|r-r'|;fc F ) V TF (r'; [p])dr' . 



(32) 



This expression is directly useful for extended systems 
where a reasonable p can be defined. For isolated sys- 
tems such as atoms and molecules where the density van- 
ishes far from the nuclei, it seems natural to replace kp 
in Eq. ([32]) by the local Fermi wavevector (LFWV) 



Mr) = (3vr 2 p(r)) 1/3 , 



(33) 



though errors may be introduced for rapidly varying den- 
sity distributions. Using Eq. ([13]) . this yields the general 
form of our proposed HQ kinetic potential: 



V H q{v- [p],k F (v)) = Vtf{t; \p\) + V w (r; [p]) 



4a 



J /(|r-r'|;fc F (r)) p 2 / 3 (r')dr'. 



(34) 



Note that the last term in Eq. 
puted in Fourier space as 



is most easily com- 



5C F 
3(2tt)= 



/(fc/2fc F (r))p 2 / 3 (k)e- 



ikr 



dk. 



(35) 



As can been seen in Eq. (|31|) . the weight function 
f(k/2k F (r)) is determined analytically. Unlike the LR- 
based KEDF approaches, no first-order differential equa- 
tion is needed to solve for the weight function in Fourier 
space. For atomic systems where the LFWV is used, the 
convolution in Eq. (|35j) must be carried out numerically, 
which will lead to a quadratic scaling of the HQ model 
(and related LQ model described below) in the number N 
of electrons. For extended systems where one can expand 
about the local density po, one can use fast Fourier trans- 
forms (FFT's) for a much more efficient computation of 
this integral (4i |. 



D. LQ Kinetic Potential 

In analogy to the TFAW model in Eq. ([2l)]) . we could 
similarly generate a KP that is accurate to second order 
at low q while still reproducing the leading term at high 
q. However this is numerically less useful because the 
analogue of Eq. (|3 5 j) involves the Fourier transform of 
Vw(r'', [p]), which cannot be simply expressed in terms 
of the density. Instead, by empirically taking a properly 
chosen component of the density outside the integral we 
find that Eq. ([34]) can be modified to produce a new LQ 
KP that is accurate to second order at low q and first 
order at high q (45|: 

V LQ (r; [p],k F (r)) = Vtf (r; \p]) + V w (r; \p\) 

+Mc F pV6( r ) //(| r -r'|;Mr)) p 1/2 (r')dr' 



(36) 

Extending these ideas we have constructed a modified 
KP that satisfies LR everywhere and is accurate to sec- 
ond order at both low and high q 45]. However the 
functional form is much more complicated, and little ad- 
ditional accuracy is gained from the improved behavior 
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at very small or very large wavevectors, since all forms 
use LR to interpolate for intermediate wavevectors, and 
this is where most errors arise in practice. Thus we will 
report results here only for the HQ and LQ models. 

V. RESULTS FOR ATOMS USING THE HQ 
AND LQ MODELS 

For completeness and to compare to earlier work, we 
first briefly discuss all-electron calculations using the pro- 
posed HQ and LQ model KPs and the full atomic poten- 
tials. We then describe results using AILPS methods. 
These are compared with the KS-DFT, the TFAW mod- 
els, and the CAT model introduced by Chacon, Alvarel- 
los, and Tarazona [13, 0, 03. The CAT model is a 
LR-based KEDF method, which gives some indication of 
shell structure. We employed the latest version, which 
uses a nonlocal two-body Fermi wave vector with a pre- 
scribed functional form depending on an empirical pa- 
rameter (3 = 1/2 (defined in Eq. (3) of Ref. [14]). This 
caused the numerical calculations [46[ to be consider- 
ably more costly than those of the LQ or HQ models, 
which used the local Fermi wave vector as in Eq. (J33J). 
All calculations are spin-restricted and use the local den- 
sity approximation (LDA) [13, IH, |49| for the exchange- 
correlation functional. 



A. All-electron calculations 

All-electron calculations consider the density response 
to the large and rapidly varying nuclear potential. Since 
the system is far from the linear response regime, quanti- 
tative results from the HQ and LQ models (or from LR- 
based KEDF methods) cannot be expected. However, by 
incorporating exact results for very large and very small 
wavevectors, these models do correct major deficiencies 
of the purely local TF model (which, e.g., predicts an 
infinite density at the nucleus!) and even give some qual- 
itative indications of atomic shell structure. 

The numerical method use the Pauli kinetic potential 
V p (r; [p]) [HHII, defined as 

V p (r;[p]) = V Ts (r;[p))-V w (r;[p]). (37) 

Since VV(r; [p]) is the exact KP for a system where the 
density can be accurately described by a single orbital, 
if V p (r; [p]) is omitted, one would essentially obtain the 
ground state density of the corresponding Boson system, 
where all the electrons are in the same orbital. If wc 
represent the full density by a single orbital function 
so that p(r) and PV(r; [p]) can be written in the forms of 
Eq. (fT5|) and Eq. (fTj?)) respectively, we can then combine 
V p (r; [p]) with the one-body potential V e ff(r; [p]) in Eq. 
(|3"|), and derive a Schrodinger-like equation for the Bose 
orbital ip(r), 

{-±V 2 + V eff (r; [p]) + V p {v- [p])} ^(r) = ^(r). (38) 



Table I: Atomic energy E using the KS, LQ, HQ, CAT and 
the TFAW models in all-electron calculations. MAPE, the 
mean absolute percentage error (relative to the KS method) of 
various OF models are given at the bottom of their respective 
columns. 





KS 


LQ 


HQ 


CAT 


TF1/5W 


TFW 


He 


-2.834 


-2.565 


-2.437 


-2.675 


-2.911 


-1.559 


Ne 


— 128.2 


—134.3 


— 126.6 


— 126.2 


—129.5 


—86.40 


Ar 


-525.9 


-545.9 


-512.2 


-515.1 


-526.2 


-375.5 


Kr 


-2750 


-2805 


-2621 


-2712 


-2748 


-2099 


Xe 


-7229 


-7306 


-6844 


-7141 


-7214 


-5701 


Be 


-14.45 


-14.39 


-13.64 


-14.11 


-14.71 


-8.699 


Mg 


-199.1 


-207.9 


-195.7 


-195.2 


-200.0 


-136.4 


C 


-37.42 


-38.97 


-36.85 


-37.25 


-38.41 


-24.01 


N 


-54.02 


-56.71 


-53.59 


-53.84 


-55.39 


-35.33 





-74.47 


-78.39 


-74.02 


-74.08 


-76.11 


-49.39 


Si 


-288.2 


-300.4 


-282.5 


-282.2 


-288.9 


-200.5 


p 


-339.9 


-354.1 


-332.7 


-332.8 


-340.6 


-238.3 


s 


-396.7 


-412.8 


-387.7 


-388.5 


-397.3 


-279.9 


MAPE 




4.06% 


3.42% 


1.82% 


1.10% 


32.0% 



Table II: Electron density at the nucleus p(0), using the KS, 
LQ, HQ, CAT and the TFAW models in all-electron calcula- 
tions. MAPE, the mean absolute percentage error (relative to 
the KS method) of various OF models are given at the bottom 
of their respective columns. 





KS 


LQ 


HQ 


CAT 


TF1/5W 


TFW 


He 


3.525 


3.088 


2.742 


3.600 


18.23 


0.9515 


Ne 


614.5 


576.6 


517.6 


613.2 


2596 


169.6 


Ar 


3819 


3642 


3282 


3812 


1.548 x 10 4 


1093 


Be 


34.86 


30.49 


27.17 


33.75 


158.2 


8.952 


Mg 


1086 


1024 


920.9 


1083 


4519 


303.0 


C 


126.0 


113.3 


101.2 


122.8 


547.9 


33.07 


N 


203.9 


185.6 


166.1 


200.0 


876.7 


54.24 


O 


308.6 


284.1 


254.6 


304.6 


1317 


83.19 


Si 


1754 


1662 


1495 


1749 


7218 


493.9 


p 


2173 


2062 


1857 


2167 


8901 


614.5 


s 


2654 


2523 


2272 


2647 


1.083 x 10 4 


753.5 


MAPE 




7.61% 


17.2% 


1.14% 


331% 


72.6% 



In other words, V e //(r; [p]) + V p (r; [p]) is now the one- 
body effective potential for the corresponding Boson sys- 
tem with the same electron density. This reduction of 
an A^-fermion problem to a Boson form is widely imple- 
mented in OF-DFT due to its numerical stability and its 
easy implementation using existing KS-DFT codes [loj . 

The associated Pauli potentials for the HQ and LQ 
models are immediately obtained by subtraction of the 
W KP from Eq. ([M]) and Eq. ([35]) respectively. The stan- 
dard finite difference method for solving Euler equations 
for the TFW models [13] are implemented for the LQ 
and HQ models, and the nonlocal terms are evaluated by 
Fourier transforms. The choices of radial grids for both of 
the real and Fourier space and other detailed numerical 
methods are given in Ref. Q • The kinetic energy for HQ 
and LQ models is computed using the Herring pathway 
inEq. ©. 

As shown in Table HI the atomic energy calculated by 
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Figure 3: (Color online) Radial density r 2 p(r) of the Kr atom 
using the KS method, the LQ and HQ models, the CAT 
model, and the TFAW models (see the inset) with the full 
nuclear potential. 
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the energy-optimized TF1/5W model is very close to the 
KS-DFT, and outperforms all the LR-based models, and 
other TFAW models. In Table [TTJ we compare the elec- 
tron density at the nucleus p(0) for various models. The 
TF1/5W model overestimates p(0) by about a factor of 4, 
while the TFW model underestimates it by about 30%. 
The predicted values of p(0) for all the LR-based mod- 
els are very close to the KS results, and are much better 
than the TFAW models. 

In Fig. ([3]), we compare the radial density distribu- 
tion r 2 p(r) of the LQ and HQ models to that predicted 
by other theories for the Kr atom. Both the TF1/5W 
and TFW models predict smooth and structureless ra- 
dial density profiles. Using the full Coulombic potential, 
all the LR-based models can predict an incipient shell 
structure for heavy atoms (Z > 30), and these results 
are typical. Since the potential is certainly far beyond 
the LR regime, these qualitative results with some sug- 
gestion of shell structure are about as good as could be 
hoped for. The surprisingly good total energies given in 
TableUfor the TF1/5W model and the LR-based models 
shows that averaged thermodynamic properties are less 
sensitive to errors in the KP than is the density profile. 
The difference in the results for the LQ, HQ and CAT 
models indicates that the LR-based OF theory is being 
used outside its range of validity. As shown below, we 
gain a significant improvement by using the AILPS to 
deal with these difficulties. 



B. Ab initio local pseudopotential calculations 

As discussed earlier, the use of pseudopotentials in non- 
local LR-based OF-DFT can improve the accuracy of the 
theory because the weaker pseudopotential is more nearly 
in the LR regime, where the theory is designed to be ac- 



Figure 4: The smooth target density p v (r) from Eq. (|40[) . with 
parameters given in Table HV1 for the Si pseudoatom used in 
the inverse-KS process, and the valence density p v (r) pre- 
dicted by the LQ and HQ models using the V pa (r) (see Fig. 
[3} corresponding to p~ v (r). The arrow indicates the location 
of r c . 

0.07 | 1 1 1 1 1 1 
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Figure 5: The AILPS V ps (r) for Si generated by the target 
density p~v(r) in Fig. [4] 
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curate. Our proposed HQ and LQ models can be used 
with any existing AILPS. However, since we want to as- 
sess the performance of these models for a wide class of 
atomic systems, we describe here a new method for de- 
termining reasonable AILPS for general atomic systems. 
These pseudopotentials will be used in all our calculations 
and can be transfered to other molecular and solid state 
environments, but we expect (and find in cases where 
comparison can be made) little change if other reason- 
able AILPS are used. 

Because of the one-to-one mapping between the effec- 
tive one-body potential acting on a system of N electrons 
and the electron density in the ground-state configura- 
tion, it is possible to obtain a unique local one-body po- 
tential that generates a given target density p(r) by using 
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Figure 6: (Color online) Radial valence density r 2 p v (r) of 
the Si atom using the KS method and various models using 
AILPS. Parameters used for constructing this reference sys- 
tem are shown in Table HV1 The arrow indicates the location 
of r c . Inset: The corresponding radial total density r 2 p(r), 
which is dominated by the core component for r < r c . 




a KS orbital-based method in an inverse way [23[ . To con- 
struct an AILPS for a given atom we separate the total 
electron density p(r) into a "core density" /5 c (r), which is 
supposed not to vary significantly in other molecular or 
solid state environments, and the target "valence density" 
p„(r) where 

p(r) = p v (v) + p c (r). (39) 

Because DFT requires only the electron density, we can 
take a more general view of what is meant by the core and 
valence components than is used in most orbital-based 
methods. Here, we directly construct a smooth target 
valence density for the N v = N — N c valence electrons, 
with N c chosen to be the number of electrons in the noble 
gas configuration. 

Our proposed target valence density p v (r) for atoms 
equals the full KS density pks{ t ) outside a core of radius 
r c , and is designed to be small and slowly varying inside 
r c . The functional form we take is 




tpKs(r c ) + a r q exp[-r p (a 1 + a 2 r 2 )} r < r c 



(40) 

Figure 2] gives an example of p v (r) for Si that will be 
discussed in more detail below. We find most results are 
insensitive to the details of our fitting procedure. Param- 
eter values for a variety of atomic systems are given in 
the Appendix in Table IIVI along with the physical and 
technical considerations that guided our choice of this 
particular form for p v (r). The Appendix also discusses 
some of the general issues that arise in using these atomic 
AILPS in other environments. 

The local pseudopotential is directly related to the 
effective one body potential that reproduces /5„(r) ex- 



actly when using the full KS theory. Following previous 
work [H,[2^], for a given p v (r), the inverse-KS equations 
are solved to get the effective one-body screened poten- 
tial V scr (r). The desired ab initio local pseudopotential 
V^ s (r) is then obtained by subtracting the Hartree po- 
tential and the exchange-correlation potential: 

V ps (r) = V scr {r) - Vb(t\ [p v ]) - V xc (r; [p v ]). (41) 

This relatively expensive procedure to determine V^ s (r) 
requires the use of orbitals. However it needs to be done 
only once for each atom, and the resulting V ps (r) can 
then be used in a variety of other environments if the 
atomic core densities remain essentially constant. 

Once suitable V ps (r) have been determined (by this or 
other means) they can be incorporated in V ex t(v) in dif- 
ferent ways, depending on the particular system of inter- 
est. OF-DFT theory can then be used to determine the 
valence density p v (r) in direct analogy to the all-electron 
calculations for the full atomic potential in Eq. (|55|) : 

{-\y 2 + V H {v- [p v ]) + V xc (r; [ Pv ]) 

+V ext (r) + V P {r; [p v ])}M r ) = WM r ) (42) 

where 

p v (r) = |^(r)| 2 . (43) 

A simple and direct test of OF-DFT is to use Eqs. 1(12)1 
and (|43]) for the same atomic system for which V ps (r) was 
constructed. Thus we take V ex t(v) = V ps (r) for a given 
atom as input data in Eq. (|4"2"]) . The valence density 
p v (r) predicted by the OF KPs is determined from Eqs. 
(|42|) and (|43|) . and can be directly compared to the exact 
target density p v (r) for this atomic system given by the 
full KS theory. 

This is illustrated in Fig. 0] which shows the input 
target valence density p v (r) for Si used in the inverse- 
KS process. The total KS density pKs(r) equals p v (r) 
for r > r c , indicated by the arrow in Fig. 21 and then 
increases rapidly for r < r c , reaching a very large value 
at the nucleus, pks(0) = 1754. In contrast, the proposed 
target valence density p v (r) remains small and relatively 
slowly varying inside the core, with p v (0) = 0.005621. 

Also shown in Fig.|4]are the predicted valence densities 
p v (r) for Si given by the LQ and HQ models. Because of 
the relatively weak V ps (r) and slowly varying valence den- 
sity, both the LQ and HQ KP models predict results very 
close to those given by the exact KS treatment of the ki- 
netic energy, and perform markedly better than they did 
for the all-electron calculations using the full Coulomb 
potential. Fig. [5] shows the corresponding AILPS gener- 
ated by the inverse KS procedure. It is much smaller in 
the core region than the full atomic potential and more 
likely to be accurately treated by LR-based methods. 

The radial valence density r 2 p v (r) of the Si and Ar 
atoms predicted by the various methods are shown in 
Figs. [5] and [71 The consistency of our OF theory when 
pseudopotentials are used is illustrated by the similar- 
ity of the density predicted by the HQ and LQ models. 
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Figure 7: (Color online) Same as in Fig. [5] but for the Ar 
atom. 




The slight deviations from the KS-DFT results for Ar in 
the valence region in Fig. [7] are among the largest we en- 
countered for all atoms tested, and could be due to the 
relatively large number of valence electrons (8) compared 
to core electrons (10). Further errors may arise from the 
LFWV approximation in Eq. . 

Once p v (r) has been determined using OF-DFT, it can 
be added to the known input core density /5 c (r) to obtain 
the predicted total density, since the basic assumption of 
our AILPS is that the core density remains unchanged 
in different chemical environments. The core density is 
defined in Eqs. ||35J) and (gUJ) . The inset in Fig. © shows 
that the HQ and LQ treatment of the valence density 
for the Si atom does not produce noticeable errors in the 
total density, as expected, and the shell structure remains 
in excellent agreement with the full KS-DFT calculations. 

As a further test of OF-DFT we can compare the en- 
ergy for the valence density given by the various methods 
to the exact valence energy for the target valence density 
determined by the inverse KS method. Table IIIII gives 
the valence energy values for the KS-DFT, and the LQ 
and the HQ models. As can been seen, both the LQ 
and HQ models give very good agreement with KS-DFT 
and perform significantly better than the other models. 
These results show that for this class of relatively weak 
pseudopotentials the OF treatment of the KP is quite 
satisfactory. 

To test of the transferability of the present AILPS, we 
also performed calculations for positive ions. Fig.[8]shows 
that the ionization energies of various atoms calculated 
using KS-DFT and the full atomic densities and those 
from the KS-DFT using the valence densities with AILPS 
are very similar. Therefore, the present AILPS are quite 
transferable to these positive ions. 

The ionization energy for models using the AILPS is 
obtained by subtracting the valence energies for systems 
with N v and N v — 1 electrons. Since the core electrons 
are assumed to be unaltered in different chemical envi- 



Table III: The total valence energy E v [p v ] using the KS 
method, the LQ and HQ models, the CAT model, and the 
TFAW models. MAE, the mean absolute error (relative to 
the KS method) of various OF models are given at the bot- 
tom of their respective columns. Parameters used in Eq. (|40p 
for such systems are given in Table ITVl 

KS LQ HQ CAT TF1/5W TFW 
Be -0.9914 -0.8955 -0.8950 -0.9583 -1.214 -0.7786 
C -6.134 -6.080 -6.100 -6.345 -7.761 -5.266 
N -11.04 -11.06 -11.09 -11.32 -13.93 -9.462 
O -18.01 -18.09 -18.09 -18.19 -22.44 -15.30" 
Si -3.771 -3.738 -3.750 -3.869 -4.467 -3.350 
P -6.474 -6.432 -6.455 -6.582 -7.385 -5.756 
S -10.20 -10.10 -10.14 -10.24 -11.27 -9.023 
Ar -21.37 -20.84 -20.91 -20.85 -22.40 -18.56" 
"MAE 0TT9 0T03 0T84 LoTO 1.312 



ronments, the atomic core energy is a constant that can- 
cels here or in other similar applications to molecules and 
solids. Limitations of the LQ and HQ models are more 
evident here, but they do capture the overall periodicity 
of the ionization energies well, and perform significantly 
better than the other models. 



VI. CONCLUSION 

In summary, we propose two nonlocal OF KPs that 
satisfy exact limits for small and large wavevector per- 
turbations and reproduce the exact LR function in the 
homogeneous limit. These are the same limits that sev- 
eral current KEDFs are designed to satisfy. However, 
because of the more local nature of the KP, it is much 
easier to satisfy these conditions for the KP than for the 
KEDF, and there may be other physical and technical 
advantages arising from use the more local KP. 

In general, there is no reason to believe that any LR- 
based OF-DFT should work well for arbitrary systems 
where the model potentials are far beyond the LR regime. 
However, most chemical processes involve changes of va- 
lence electron densities, which can often be described by 
a weak AILPS. Thus the use of a LR-based OF-KP to- 
gether with AILPS for atomic systems at least seems well 
justified. The small and relatively slowly varying p v (r) 
also provides some justification for our use of the local 
FWV k F {r) in Eq. 

When the AILPS is used, the valence densities given by 
the LQ and HQ KPs are close to those given by the KS 
method. Thus the particular integration pathway used 
to get the total energy value becomes unimportant. The 
simple pathway in Eq. 1(5)) is especially useful, since no 
coupling parameter integration is needed. 

The proposed models are not only conceptually simple, 
but also exact for a model system with a weak poten- 
tial and a slowly-varying density. The appearance of the 
atomic shell structure was found to be very sensitive to 
the accuracy of the proposed KPs. The LR-based LQ and 



Figure 8: (Color online) Ionization energies (shown in 
Hartree) of the first and the second row atoms using the full 
KS method, and various models using AILPS. The mean abso- 
lute errors (relative to the full KS method) of various models 
using AILPS are: KS (0.1 eV), LQ (1.8 eV), HQ (2.1 eV), 
CAT (3.9 eV), TFW (3.1 eV), and TF1/5W (4.3 eV). Ioniza- 
tion energies using the TF1/5W model are not shown in the 
figure due to its relatively poor performance. 



0.9 




0.1 1 ' 1 1 ' ' 1 1 ' ' 1 1 ' ' 1 1 

3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 

z 



HQ KPs give at best only qualitative indications of shell 
structure for full atomic systems, though total energies 
are surprisingly good. Still better results for atoms and 
ions can be found by focusing on the valence density as 
determined by a relatively weak AILPS. While these re- 
sults seem promising, improved KPs are needed and fur- 
ther investigation is required to see if these ideas can be 
usefully applied to other relevant systems like molecules 
and solids. Some initial results along these lines will be 
reported elsewhere. 
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Table IV: Parameters used in Eq. (|A.1[) for the target valence 
density of various atoms. Here, p — q — 6, and t = 0.1 
are used for all systems. The V pa (r) generated from these 
parameterized p v (r) can then be used in OF-DFT. 
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Ne 


8 


1.052 x 10 


4 


1993 


-8267 


0.3707 


Na 


1 


5.234 x 10~ 


b 9. 


840 x 10" d 


-6.610 x 10~ 4 


2.904 


Mg 


2 


9.805 x 10~ 


4 


0.04445 


-5.002 x 10"^ 


2.233 


Al 


3 


6.164 x 10" 


y 


0.1273 


-0.02056 


1.861 


Si 


4 


0.02942 




0.3119 


-0.06861 


1.593 


P 


5 


0.1149 




0.6864 


-0.1979 


1.390 


S 


6 


0.3841 




1.390 


-0.5100 


1.231 


CI 


7 


1.133 




2.632 


-1.201 


1.103 


Ar 


8 


3.018 




4.714 


-2.624 


0.9985 



Appendix: AB INITIO LOCAL 
PSEUDOPOTENTIALS 

As discussed above, our proposed target valence den- 
sity for atoms has the following form: 




tpKs(r c ) + a r q exp[-r p (a 1 + a 2 r 2 )} r < r c 



Here p and q are taken as even integers. The larger 
they are, the smaller and more slowly varying is the va- 
lence density near r — but the sharper is the peak 
near the core radius r c . As a compromise, we take here 
p = q = 6, which generates relatively slowly varying local 
pseudopotentials V ps (r). For applications in different en- 
vironments, such as molecules or crystals, the core size r c 
has to be small to maintain transferability of the atomic 
core density. For this reason, we force p v [r — 0) to be 
small by taking a small t. If t = 0, the strict vanish- 
ing of the valence density near the nucleus would require 
a very repulsive V ps , which is certainly undesirable for 
the LR-based OF-DFT. However, if t is too large, there 
will exist a long oscillatory tail outside the core in the 
corresponding V ps (r). This is an undesirable feature for 
transferability to other environments, as will be discussed 
below. These two points constrain the value of t and we 
use here t = 0.1 for all the atomic systems considered. 

The four parameters eto, a\, 0,2, and r c are determined 
by requiring continuity of the function p v (r) and its first 
two derivatives at r = r c , and by satisfying the normal- 
ization condition: 

N v = 4tt y p v (r)r 2 dr. (A.2) 
Here we used the standard noble gas cores to determine 
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Figure 9: The inverse KS procedure generates very small 
oscillations in the tail of rV p3 (r) for Si (shown here) and 
other atoms. The two points where rV pa {r) — —4 for Si are 
n = 2.336 and r 2 = 4.576 or r\ /2 = 2.139. The arrow indi- 
cates the location of r 2 . To achieve good transferability (see 
text), this rV ps (r) is modified by setting rV ps (r) = —4 when 
r > T2 = 4.576. See Fig. [5] for a large scale view. 
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N v = N — N c , though other choices could in principle be 
made. 

To construct our local pseudopotential V ps (r) for 
atoms, we first solve the KS equations for an atom with 
the full Coulomb potential. With the KS density pKs{ r ) 
and N c determined, this construction ensures that as 
r — > r c , the associated core density smoothly approaches 
zero as 0(| r — r c | 3 ). 

After generating the parameterized target valence den- 
sity p v (r), the set of inverse-KS equations are solved 
to obtain the corresponding one-body screened poten- 
tial V scr (r). In principle, the AILPS V ps (r) is then given 
by Eq. (|4"Tj) . Using this, we find with an acceptable N c 
that essentially the same density profiles outside r c are 
predicted by the LQ and HQ models for a wide range of 
choices of p, q, and t. This would be expected if most fea- 
tures of the resulting set of model potentials are within 
the different regimes accurately described by the OF KP. 



These results clearly show that OF-DFT can give accu- 
rate results for this class of physically relevant and rela- 
tively weak and slowly varying model potentials. 

However, we found that the rV ps {r) constructed in this 
way for group III to group VIII elements deviates from 
— Z v outside the core by a small and long range oscilla- 
tion. For the 2s 2 2p x atoms, the maximum amplitude of 
the oscillation is about 0.18. For the 3s 2 3p x atoms, it 
ranges from 0.003 for Al to 0.037 for Ar. Similar oscil- 
latory tails were observed by Wang and Stott [25J, and 
are thought to arise from the inability of a local pseu- 
dopotential to represent both s and p orbitals of the cor- 
responding nonlocal pseudopotentials. These oscillations 
are so small in magnitude that they have almost no effect 
on the density profile or energy of atoms, but they can 
cause transferability problems when used for solids [o3l |. 
Following Wang and Stott's approach, we replace V ps (r) 
by —Z v /r at the largest point where rV ps (r) = ~Z V (see 
Fig. |9])). The new V ps (r) can then be used for both 
atomic and solid state calculations. 

In other chemical environments, such molecules and 
solids, the predetermined AILPS V ps (r) centered at each 
nucleus are regarded as input data, and the valence den- 
sities for other systems can be calculated by OF-DFT. 
Results can be checked by using the full KS-DFT. For 
example, for molecules and solids, the external potential 
Vexti?) is a linear combination of the local atomic pseu- 
dopotentials centered at each ion position R/: 

Wr)=5Z^ p .(r-R / ) (A.3) 

The valence density can then be determined from Eqs. 
(|42)l and ([43]) . and the total energy of valence electrons 
is 

K\fiv] = T s [p v ] + E H [p v ] + E xc [p v ] + J p v (r)V ext {v)dv 

(A.4) 

For the proposed LQ and HQ KPs, the value of T^fp^] 
in Eq. (|A.4p is determined by Eqs. (J5J or (fTT|) with p(r) 
replaced by p v (r). 
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